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Grid-independent  large-eddy  simulation  in 
turbulent  channel  flow  using  three-dimensional 

explicit  filtering 

By  Jessica  Gullbrand 


1.  Motivation 

The  most  commonly  used  Large  Eddy  Simulation  (LES)  approach  is  the  implicitly 
filtered  approach.  In  implicitly  filtered  LES,  the  computational  grid  and  the  discretization 
operators  are  considered  as  the  filtering  of  the  governing  equations.Thereby  the  turbulent 
flow  field  is  divided  into  grid  resolved  and  unresolved  scales,  where  the  unresolved  scales 
must  be  modeled. 

When  explicit  filtering  is  used  in  LES,  the  filtering  procedure  of  the  governing  equations 
is  separated  from  the  grid  and  discretization  operations.  The  flow  field  is  divided  into 
resolved  filtered  scale  (RFS)  motions,  and  subfilter-scale  (SFS)  motions.  The  SFS  is  itself 
divided  into  a resolved  part  (RSFS)  and  an  unresolved  part  (USFS)  (Zhou  et  al.  2001); 
see  figure  1.  The  RFS  motion  is  obtained  by  solving  the  filtered  Navier-Stokes  equations. 
The  RSFS  motions  can  be  reconstructed  from  the  resolved  field  and  occur  due  to  the  use 
of  a smooth  (in  spectral  space)  filter  function.  The  USFS  motions  consist  of  scales  that  are 
not  resolved  in  the  simulation  and  need  to  be  modeled.  The  explicitly  filtered  governing 
equations  were  recently  studied  by  Carati  et  al.  (2001)  in  forced  isotropic  turbulence. 

The  smallest  resolved  scales  are  often  used  to  model  the  turbulence- closure  term  in 
LES,  and  therefore  it  is  important  to  capture  these  scales  accurately.  The  accuracy  of  the 
LES  solution  can  be  increased  by  using  high-order  numerical  schemes.  Although  high- 
order  methods  treat  the  important  large  energy-containing  scales  more  accurately,  the 
small  resolved  scales  would  still  be  contaminated  with  truncation  errors  when  using  non- 
spectral methods.  These  errors  can  be  reduced  or  eliminated  by  using  explicit  filtering 
in  LES  (Lund  1997).  This  can  be  achieved  either  by  using  a large  ratio  of  filter  width  to 
cell  size,  or  by  using  a higher-order  method,  in  which  case  the  ratio  need  not  be  so  large. 
In  recent  a priori  studies  by  Chow  & Moin  (2003),  a minimum  ratio  of  filter  width  to 
cell  size  was  determined  to  prevent  the  numerical  error  from  becoming  larger  than  the 
contribution  from  the  turbulence-closure  term.  They  concluded  that  with  a fourth-order 
scheme  a filter  width  of  at  least  twice  the  cell  size  should  be  used,  and  for  a second-order 
scheme  the  filter  width  should  be  at  least  four  times  the  cell  size. 

Using  explicit  filtering  and  high-order  numerical  schemes  requires  the  filter  functions 
to  be  commutative  to  at  least  the  same  order  as  the  numerical  scheme.  The  differenti- 
ation and  filtering  operations  must  commute,  to  ensure  that  the  filtered  Navier-Stokes 
equations  have  nearly  the  same  structure  as  the  unfiltered  equations.  In  general,  the 
operations  do  not  commute  when  a variable  filter  width  is  used,  as  is  needed  in  inho- 
mogeneous turbulent  flows.  Ghosal  k Moin  (1995)  showed  that  the  commutation  error 
might  overwhelm  the  contribution  from  the  turbulence-closure  term.  Therefore,  this  error 
must  be  reduced  or  eliminated  to  avoid  significant  effects  on  the  LES  solution.  A general 
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theory  for  constructing  discrete  high-order  commutative  filters  was  proposed  by  Vasilyev 
et  al.  (1998). 

Most  of  the  previous  studies  of  LES  using  explicit  filtering  in  turbulent  channel  flow 
have  used  filtering  in  two  dimensions  (the  homogeneous  directions)  and  only  a few  studies 
have  applied  filtering  in  all  three  dimensions.  We  limit  our  discussion  to  investigations 
performed  using  smooth  filter  functions.  Two-dimensional  filtering  was  investigated  by 
Moin  k.  Kim  (1982),  Piomelli  et  al.  (1988),  Najjar  & Tafti  (1996),  and  Gullbrand  & Chow 
(2002)  among  others.  Studies  using  three-dimensional  filtering  were  performed  by  Cabot 
(1994),  Gullbrand  (2001),  Winckelmans  et  al.  (2001),  and  Stolz  et  al.  (2001).  However, 
most  of  the  studies  using  three-dimensional  filtering  did  not  focus  on  minimizing  the 
effect  of  the  numerical  errors.  If  care  is  not  taken  to  reduce  the  numerical  errors,  they 
may  be  larger  than  the  contribution  from  the  turbulence-closure  models.  Therefore,  it 
will  not  be  possible  to  separate  the  numerical  effects  from  the  performance  of  turbulence- 
closure  models.  Cabot  (1994),  for  example,  used  a second-order  finite-difference  scheme 
and  second-order  commutative  filter  functions  with  a ratio  of  two  between  the  local 
filter  width  and  the  local  cell  size.  The  error  from  the  second-order  scheme  is  proba- 
bly larger  than  the  turbulence-closure  contribution  due  to  the  small  ratio  of  the  filter 
width  to  cell  size  used  and,  in  addition,  a second-order  commutation  error  is  present. 
Winckelmans  et  al.  (2001)  used  a high-order  finite-difference  scheme  (fourth-order)  but 
applied  a second-order  commutative  filter  with  a ratio  of  filter  width  to  cell  size  of  \/6. 
The  filter  functions  used  introduce  a commutation  error  of  second-order  into  the  simu- 
lations. A spectral  method  was  used  by  Stolz  et  al.  (2001),  together  with  fourth-order 
commutative  filter  functions  with  a filter-grid  ratio  of  approximately  1.5.  The  use  of  spec- 
tral methods  clearly  reduces  the  numerical  errors  in  the  simulation  when  compared  to 
the  studies  previously  mentioned.  However,  spectral  methods  axe  not  considered  in  this 
study.  Gullbrand  (2001)  used  fourth-order  commutative  filter  functions,  with  a ratio  of 
two  between  the  local  filter  width  and  the  local  cell  size,  in  a fourth-order  finite-difference 
code.  The  commutation  error  is  then  of  the  same  order  as  the  numerical  scheme,  which 
is  of  higher  order  than  the  turbulence-closure  contribution.  According  to  the  study  by 
Chow  & Moin  (2003),  the  filter-grid  ratio  used  ensures  that  the  contribution  from  the 
turbulence-closure  term  is  larger  than  the  numerical  errors.  Thus,  a fourth-order  scheme 
using  fourth-order  commutative  filters  with  a filter  width  of  at  least  twice  the  cell  size 
creates  a numerically-clean  environment  where  turbulence-closure  models  can  be  tested 
and  validated. 

In  this  paper,  turbulence-closure  models  are  evaluated  using  the  “true”  LES  approach 
in  turbulent  channel  flow.  The  study  is  an  extension  of  the  work  presented  by  Gullbrand 
(2001),  where  fourth-order  commutative  filter  functions  are  applied  in  three  dimensions 
in  a fourth-order  finite- difference  code.  The  true  LES  solution  is  the  grid-independent 
solution  to  the  filtered  governing  equations.  The  solution  is  obtained  by  keeping  the  filter 
width  constant  while  the  computational  grid  is  refined  (figure  2).  As  the  grid  is  refined, 
the  solution  converges  towards  the  true  LES  solution.  The  true  LES  solution  will  depend 
on  the  filter  width  used,  but  will  be  independent  of  the  grid  resolution.  In  traditional 
LES,  because  the  filter  is  implicit  and  directly  connected  to  the  grid  spacing,  the  solution 
converges  towards  a direct  numerical  simulation  (DNS)  as  the  grid  is  refined,  and  not 
towards  the  solution  of  the  filtered  Navier-Stokes  equations.  The  effect  of  turbulence- 
closure  models  is  therefore  difficult  to  determine  in  traditional  LES  because,  as  the  grid 
is  refined,  more  turbulence  length  scales  are  resolved  and  less  influence  from  the  models  is 
expected.  In  contrast,  in  the  true  LES  formulation,  the  explicit  filter  eliminates  all  scales 
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that  are  smaller  than  the  filter  cutoff  {keg  in  figure  2),  regardless  of  the  grid  resolution. 
This  ensures  that  the  resolved  length-scales  do  not  vary  as  the  grid  resolution  is  changed. 
In  true  LES,  the  cell  size  must  be  smaller  than  or  equal  to  the  cutoff  length  scale  of  the 
filter  function. 

The  turbulence-closure  models  investigated  are  the  dynamic  Smagorinsky  model  (DSM), 
the  dynamic  mixed  model  (DMM),  and  the  dynamic  reconstruction  model  (DRM).  These 
turbulence  models  were  previously  studied  using  two-dimensional  explicit  filtering  in  tur- 
bulent channel  flow  by  GuUbrand  & Chow  (2002).  The  DSM  by  Germano  et  al.  (1991) 
is  used  as  the  USES  model  in  all  the  simulations.  This  enables  evaluation  of  different 
reconstruction  models  for  the  RSFS  stresses.  The  DMM  (Zang  et  al.  1993)  consists  of 
the  scale-similarity  model  (SSM)  by  Bardina  et  al.  (1983),  which  is  an  RSPS  model,  in 
linear  combination  with  the  DSM.  In  the  DRM  (GuUbrand  & Chow  2002),  the  RSFS 
stresses  are  modeled  by  using  an  estimate  of  the  unfiltered  velocity  in  the  unclosed  term 
(Stolz  et  al.  2001),  while  the  USFS  stresses  are  modeled  by  the  DSM.  The  DSM  and 
the  DMM  are  two  commonly  used  turbulence-closure  models,  while  the  DRM  is  a more 
recent  model. 


2.  Governing  equations 

The  governing  equations  for  an  incompressible  flow  field  are  the  continuity  equation 
together  with  the  Navier-Stokes  equations, 

^=0  guj  dujUj  ^ dp  1 d'^Uj 

dxi  ’ dt  dxj  dxi  Rct  dxjdxj  ’ ^ ’ 

Here  Ui  denotes  velocity,  p pressure  and  Rcr  the  Reynolds  number  based  upon  friction 
velocity,  Ur,  and  channel  half-width,  h.  Einstein  summation  is  applied  to  repeated  indices. 

In  LES,  the  governing  equations  are  filtered  in  space.  The  filtering  procedure  is  applied 
to  the  flow-field  variables  according  to 

Ui{x,A,t)=  I G{x,x',A)ui{x',t)dx' , (2.2) 

Jd 

where  G is  the  filter  function  and  A is  the  filter  width. 

Hence,  the  filtered  governing  equations  can  be  written  as 


dui  _ Q dui  duiUj  _ dp  1 d'^Ui  dfij 

dxi  ' dt  ^ dxj  dxi  ^ Rct  dxjdxj  dxj 


(2.3) 


where  the  turbulent  stresses  are  defined  as  Ty  = ujuj  — %Uj.  The  filtered  equations  are 
not  closed  because  of  the  nonlinear  term  uiuj.  The  approach  by  Leonard  is  followed  where 
the  instantaneous  velocity  field  is  divided  into  a filtered  velocity  and  a fluctuating  part, 
Ui  = Ui-hu'i,  and  the  unclosed  term  can  be  rewritten  to  = (uj  -f  u'i){uj  -f  u' ) = 
UiUj  d-  Tij.  The  product  of  the  nonlinear  terms  {uiUj)  introduces  high  wavenumbers 
that  are  beyond  the  wavenumber  content  of  the  filtered  velocity  field  (uj)  and  beyond 
the  cutoff  of  the  filter  function.  To  prevent  these  high  wavenumbers  to  influence  the 
resolved  wavenumbers,  the  nonlinear  terms  are  explicitly  filtered.  A potential  drawback 
of  (2.3)  is  that  the  resulting  equation  is  not  in  general  Galilean-invariant  (Speziale  1985) 
provided  that  the  SFS  models  are  Galilean-invariant.  In  a moving  coordinate  system,  the 
resulting  equation  contains  the  additional  term  Cjd^i  — Ui)ldxj,  where  Cj  is  the  uniform 
translation  velocity.  The  additional  term  is  proportional  to  the  difference  between  the 
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doubly  filtered  and  singly  filtered  velocity  field.  This  difference  will  be  zero  when  a sharp 
cut-off  filter  is  used,  but  will  not  vanish  in  the  general  case.  The  error  can  be  minimized 
by  constructing  the  explicit  filter  function  as  close  as  possible  to  a sharp  cut-off  filter. 
However,  the  preferable  solution  is  to  choose  an  appropriate  turbulence-closure  model  so 
the  problem  can  be  avoided.  Speziale  (1985)  showed  that  using  the  SSM  by  Bardina  et  al. 
(1983)  with  the  model  coefficient  of  unity  as  the  turbulence- closure  model  would  solve  the 
problem.  However,  a drawback  of  this  model  is  that  it  generates  higher  frequencies  and 
the  desired  effects  from  the  explicit  filtering  approach  is  thus  destroyed  (Lund  1997).  The 
DRM  on  the  other  hand  has  the  desired  properties  and  avoids  the  Galilean-invariance 
problem.  Description  of  this  model  and  further  discussion  are  found  in  the  next  section. 


3.  Subfilter-scale  models 

The  turbulent  flow  field  is  divided  into  RFS  and  SFS  motions  when  explicit  filtering  of 
the  Navier-Stokes  equations  is  applied.  In  figure  1,  a sketch  of  a typical  energy  spectrum 
is  shown.  The  solid  line  represents  the  energy  captured  by  a fully  resolved  DNS,  while 
the  dashed  line  represents  the  LES  energy.  The  vertical  line  at  keg  shows  the  filter  cutoff 
in  the  LES.  The  filter  cutoff  is  determined  by  where  the  filter  function  goes  to  zero  and 
stays  zero  (in  spectral  space),  i.e.,  no  wavenumbers  higher  than  the  cutoff  wavenumber 
are  resolved  in  the  simulation.  The  filter  cutoff  can  be  seen  in  figure  2.  All  wavenumbers 
smaller  than  the  filter  cutoff  wavenumber  are  resolved  in  the  simulations.  However,  they 
are  damped  by  the  filter  function  and  have  to  be  recovered  by  an  inverse  filter  operation. 
This  corresponds  to  the  RSFS  portion  of  the  energy  spectrum.  The  same  terminology  for 
the  RSFS  and  the  USFS  was  previously  suggested  by  Zhou  et  al.  (2001).  In  principle,  the 
RSFS  can  be  exactly  recovered,  but  this  is  only  possible  when  using  spectral  methods. 
If  non-spectral  methods  are  applied,  there  are  numerical  errors  (NE)  associated  with  the 
high  wavenumbers  and  thus  the  recovered  scales  are  contaminated  with  errors. 

The  unresolved  portion  of  the  spectrum  (the  USFS)  consists  of  wavenumbers  that  are 
higher  than  the  filter  cutoff  wavenumber.  The  USFS  motions  need  to  be  modeled.  The 
vertical  lines  in  figure  1 represent  the  grid  cutoff  wavenumbers  for  two  grid  resolutions. 
The  coarse  grid  cutoff,  keg,  happens  to  coincide  with  the  filter  cutoff,  while  the  fine  grid 
cutoff,  kfg,  is  located  in  the  USFS  portion  of  the  spectrum.  However,  the  USFS  motions 
that  need  to  be  modeled  are  the  same  for  the  two  resolutions,  since  the  filter  cutoff 
determines  the  wavenumbers  resolved. 

To  recover  the  RSFS  stresses,  the  iterative  method  of  van  Cittert  (1931)  is  used  in 
this  study.  This  method  was  previously  used  by  Stolz  et  al.  (2001)  in  their  approximate 
deconvolution  procedure  to  reconstruct  the  unfiltered  velocity  field  Uj  from  the  filtered 
field  Ui.  To  fully  recover  the  unfiltered  velocity,  an  infinite  number  of  iterations  is  needed. 
However,  since  this  is  not  practical  in  numerical  simulations,  the  unfiltered  velocity  field 
is  approximated  by  a finite  number  of  iterations.  By  varying  this  number,  different  models 
can  be  obtained  to  model  the  RSFS  stresses. 

Here,  low-level  reconstruction  (the  SSM)  and  reconstruction  up  to  level  five  are  used. 
Further  details  of  the  reconstruction  used  are  found  in  3.3.  In  order  to  compare  the 
different  RSFS  models,  the  same  USFS  model  (the  DSM)  is  used  in  all  the  simulations. 
The  combinations  of  RSFS  and  USFS  models  used  are  described  below. 
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Figure  1.  Schematic  of  velocity  energy  spectrum  showing  partitioning  into  resolved  filtered 
scale  (RFS),  resolved  subfilter-scale  (RSFS),  and  unresolved  subfilter-scale  (USFS)  motions.  The 

numerical  error  (NE)  region,  denoted  by , is  a subregion  of  the  RSFS. represents 

DNS  energy, LES  energy,  and filter  cutoff.  The  vertical  line  at  keg  represents  the 

filter  cutoff  wavenumber,  which  corresponds  to  the  smallest  resolved  wavenumber  for  the  coarse 
grid.  The  vertical  line  at  k/g  represents  the  wavenumber  cutoff  for  the  fine  grid. 

3.1.  Dynamic  Smagorinsky  Model 

The  DSM  is  a widely-used  eddy  viscosity  USFS  model  (Smagorinsky  1963): 

Tig  = -2i^ij  = -2(CA)2|S|5,^.,  (3.1) 

where  i/g  is  the  eddy  viscosity,  A the  filter  width,  and  the  strain  rate  tensor  = 
0.5{dui/dxj  + duj/dxi).  The  model  parameter  (CA)^  is  calculated  dynamically  (Ger- 
mane et  al.  1991)  using  the  least-square  approximation  of  Lilly  (1992).  The  model  pa- 
rameter is  averaged  in  the  homogeneous  directions  and  is  calculated  by  the  same  dynamic 
procedure  as  described  in  the  papers  previously  mentioned.  The  explicit  filtering  of  the 
nonlinear  terms  is  not  considered  when  the  model  parameter  is  calculated.  The  filter- 
ing enters  only  when  is  introduced  into  the  filtered  Navier-Stokes  equations.  Large 
magnitudes  of  negative  values  of  the  eddy  viscosity  are  clipped  to  avoid  negative  total 
viscosity  in  the  simulations,  {i/e  + i'  > 0),  as  proposed  by  Zang  et  al.  (1993). 

3.2.  Dynamic  Mixed  Model 

Low-level  reconstruction  of  the  RSFS  stresses  can  be  performed  by  using  the  SSM  pro- 
posed by  Bardina  et  al.  (1983).  Here  the  RSFS  stress  is  modeled  by  the  scale-similarity 
term  and  the  DSM  is  used  as  the  USFS  model: 

Tig  = (S^J-  - iilj)  - 2(C'A)21SlSij.  , (3.2) 

to  form  the  DMM.  The  SSM  term  is  discretized  with  the  same  numerical  scheme  as  the 
convective  terms. 


3.3.  Dynamic  Reconstruction  Model 

High-order  reconstruction  of  the  RSFS  stress  tensor  can  be  achieved  by  the  iterative 
deconvolution  method  of  van  Cittert  (1931).  The  unfiltered  quantities  can  be  derived  by 
a series  of  successive  filtering  operations  (G)  applied  to  the  filtered  quantities  with 

Ui  —Ui  + {I  - G)  *Ui  -{■  (I  - G)  * ((/  -G)*Ui)-] (3.3) 

where  I is  the  identity  matrix.  The  truncation  order  of  the  expansion  determines  the 
level  of  deconvolution,  as  discussed  by  Stolz  et  al.  (2001).  If  the  series  includes  the  terms 
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explicitly  shown  in  (3.3),  it  corresponds  to  reconstruction  of  level  two.  An  approximate 
unfiltered  velocity  (u*)  is  obtained  by  the  truncated  series,  u*  is  substituted  into  the 
unclosed  term  which  results  in  The  reconstruction  of  the  RSFS  stresses  are 
used  in  linear  combination  with  the  DSM, 

Tij  = u*^u*j  - - 2(CA)2|5'|5^j.  , (3.4) 

which  is  called  the  dynamic  reconstruction  model  (DRM).  In  the  simulations,  the  same 
numerical  scheme  is  used  for  the  convective  terms  and  the  RSFS  terms.  The  DRM  yields 
a Galilean-invariant  expression  of  (2.3),  since  the  nonlinear  terms  UiUj  on  the  right-hand 
side  and  left-hand  side  of  the  equation  cancel  each  other.  A reconstruction  series  of  up 
to  level  five  is  used  in  this  study. 


4.  Filter  functions 

It  is  important  that  the  explicit  filter  and  the  test  filter,  which  is  used  in  the  dynamic 
procedure  of  the  DSM,  have  similar  shapes,  since  the  dynamic  procedure  is  based  upon 
the  scale-similarity  assumption  in  the  Germano  identity  (Germano  et  al.  1991).  In  the 
simulations  presented  here,  the  same  filter  function  is  used  in  all  the  simulations.  It  is 
only  the  filter  width  that  is  varied  between  the  simulations.  The  base  filter  is  a fourth- 
order  commutative  filter  function  with  filter  width  2Acg,  where  Acp  is  the  grid  cell  size 
for  the  coarse-grid  resolution.  The  computational  domain  and  grid  resolutions  used  in 
the  simulations  are  discussed  in  section  6.  It  is  not  straightforward  to  determine  the 
filter  width  of  a high-order  filter,  and  different  methods  were  studied  by  Lund  (1997). 
Here,  one  of  the  methods  suggested  by  Lund  is  applied.  The  filter  width  is  defined  as  the 
location  where  the  filter  function  reaches  a value  of  G{k)  = 0.5.  The  filter  function  used 
in  the  simulations  was  developed  by  Vasilyev  et  al.  (1998)  and  is 


(4.1) 


where  the  filter  weights  for  (f)i±-2  are  zero.  The  smooth  filter  function  is  shown  in  spectral 
space  in  figure  2.  In  the  near- wall  region,  asymmetric  filters  are  used  in  the  first  three 
grid  points  for  the  coarse  grid  in  the  wall-normal  direction.  Since  the  filter  width  is  held 
fixed,  this  corresponds  to  using  asymmetric  filters  in  the  first  six  grid  points  for  the  fine 
grid. 

In  the  simulations,  the  ratio  of  the  test-filter  width  to  the  explicit-filter  width  is  chosen 
to  be  two,  as  proposed  by  Germano  et  al.  (1991)  for  the  DSM.  The  test  filter  is  used 
only  in  the  calculation  of  (CA)^  in  the  DSM,  while  the  explicit  filter  function  is  used  to 
determine  the  RSFS  contribution  through  either  the  SSM  or  reconstruction  by  the  van 
Cittert  (1931)  iterative  method.  The  ratio  between  the  explicit-filter  width  and  the  cell 
size  for  the  coarse  grid  is  two  and  for  the  fine  grid,  the  ratio  is  four.  This  preserves  the 
effective  filter  width  as  the  grid  resolution  is  increased,  as  seen  in  figure  2.  The  vertical 
line  at  low  wavenumber  represents  the  grid  cutoff  {k^g)  for  the  coarse-grid  resolution.  The 
filter  cutoff  wavenumber  is  the  same  as  the  grid  cutoff  for  the  coarse  grid.  For  the  fine 
grid,  the  filter  cutoff  is  held  fixed,  resulting  in  a separation  between  the  filter  cutoff  and 
grid  cutoff  (kfg)  locations.  The  grid-cutoff  wavenumbers  are  also  shown  schematically  in 
figure  1.  The  ratio  of  two  for  the  filter  width  and  the  cell  size  for  the  coarse  grid  was 
chosen  to  prevent  the  numerical  errors  fi:om  becoming  larger  than  the  contribution  of  the 
turbulence-closure  model  (Ghosal  1996;  Chow  & Moin  2003). 
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Figure  2.  The  base  filter  function  in  spectral  space  and  its  relation  to  the  computational  grid 

resolutions.  : filter  function,  : filter  width, : coarse  grid  resolution,  and 

: fine  grid  resolution.  The  wavenumber  keg  represents  the  grid  cutoff  wavenumber  for  the 

coarse  grid,  while  kfg  represents  the  grid  cutoff  wavenumber  for  the  fine  grid. 

5.  Solution  algorithm 

In  the  computational  code,  the  spatial  derivatives  are  discretized  using  a fourth-order 
central-difference  scheme  on  a staggered  grid.  The  convective  term  is  discretized  in  the 
skew-symmetric  form  (Morinishi,  Lund,  Vasilyev  k Moin  1998;  Vasilyev  2000)  to  ensure 
conservation  of  turbulent  kinetic  energy.  The  equations  are  integrated  in  time  ii.sing 
the  third-order  Runge-Kutta  scheme  described  by  Spalaxt,  Moser  k Rogers  (1991).  The 
diffusion  terms  in  the  wall-normal  direction  are  treated  implicitly  by  the  Crank-Nicolson 
scheme.  The  splitting  method  of  Dukowicz  k Dvinsky  (1992)  is  used  to  enforce  the 
solenoidal  condition.  The  resulting  discrete  Poisson  equation  for  the  pressure  is  solved 
in  the  wall-normal  direction  using  a penta^diagonal  matrix  solver.  In  the  homogeneous 
directions,  the  Poisson  equation  is  solved  using  a discrete  Fourier  transform.  Periodic 
boundary  conditions  are  applied  in  the  streamwise  and  spanwise  homogeneous  directions, 
and  no-slip  conditions  are  enforced  at  the  channel  walls.  A fixed  mean  pressure  gradient  is 
used  to  drive  the  flow.  The  results  using  the  fourth-order  computational  code  is  compared 
to  results  from  a second-order  finite-difference  code  in  Gullbrand  (2000)  and  Gullbrand 
k Chow  (2002). 


6.  Turbulent  channel  flow  simulations 


The  Reynolds  number  is  Re-r  = 395  and  the  computational  domain  is  (27rh,  2h,  tt/i) 
in  (x,  y,  z),  where  x is  the  streamwise  direction,  y the  wall-normal  direction,  and  z the 
spanwise  direction.  The  computational  grid  is  stretched  in  the  y-direction  by  a hyperbolic 
tangent  function 


where  N2  is  the  number  of  grid  points  in  the  wall-normal  (j)  direction  and  7 is  the 
stretching  parameter,  which  is  set  to  2.75.  Two  computational  grids  are  used  for  the 
LES  calculations;  (64,49,48),  which  corresponds  to  one-quarter  of  the  DNS  resolution  in 
each  spatial  direction,  and  (128,97,96),  which  is  half  the  number  of  DNS  grid  points  in 
each  direction.  The  cell  size  for  the  coarser  grid  resolution  is  Ax+  = 39,  = 26,  and 

0.4  < Ay+  < 45.  The  finer  resolution  corresponds  to  the  cell  size  Ax+  = 19,  Az+  = 13, 
and  0.2  < Ay+  < 23.  The  ‘plus’  values  (wall  units)  are  obtained  by  normalizing  the  length 
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Figure  3.  Mean  velocity  profiles  using  different  turbulence-closure  models,  o:  filtered  DNS, 

; DSM  (64,49,48), : DSM  (128,97,96), : DMM  (64,49,48), ; DMM 

(128,97,96), : DRM  (N=5)  (64,49,48),  and : DRM  (N=5)  (128,97,96). 


scale  with  the  friction  velocity  and  the  kinematic  viscosity.  A statistically  stationary 
solution  is  obtained  after  30  dimensionless  time  units,  and  thereafter  statistics  were 
sampled  during  15  additional  time  units.  The  time  is  normalized  with  the  friction  velocity 
and  channel  half-width.  The  LES  results  are  compared  to  filtered  DNS  data,  and  all  the 
presented  results  are  averaged  in  the  homogeneous  directions. 

The  DNS  is  performed  with  the  same  computational  code  as  used  in  the  LES  sim- 
ulations. The  computational  grid  resolution  is  (256,193,192)  as  used  by  Moser,  Kim  & 
Mansour  (1999).  The  DNS  data  fields  are  filtered  using  the  commutative  filter  functions. 
In  the  figures,  it  is  only  every  third  grid  point  in  the  filtered  DNS  that  are  plotted  to 
make  the  comparison  between  the  LES  and  DNS  clearer. 


7.  Results 

Figure  3 shows  mean  velocity  profiles  firom  simulations  using  different  RSFS  models  and 
different  grid  resolutions.  The  filter  width  is  fixed,  while  the  grid  resolution  is  increased. 
The  goal  is  to  obtain  a grid-independent  LES  solution  so  that  the  behavior  of  turbulence- 
closure  models  can  be  evaluated.  The  changes  in  the  predicted  mean  velocity  profiles  as 
the  grid  resolution  is  increased  are  only  minor,  indicating  that  the  LES  solutions  are 
nearly  converged.  The  mean  velocities  predicted  by  the  DSM  are  much  higher  than  the 
filtered  DNS  results.  The  DMM  improves  the  results  slightly,  while  the  best  agreement 
with  the  filtered  DNS  data  is  predicted  by  the  DRM  (N=5).  This  shows  the  need  for  a 
RSFS  model  when  a smooth  explicit  filter  function  is  applied. 

The  streamwise  velocity  fluctuations  in  figure  4 show  the  same  trend  as  the  mean 
velocity  profiles.  However,  the  differences  in  the  results  as  the  grid  is  refined  are  slightly 
larger  than  for  the  mean  velocity.  The  DSM  shows  the  largest  overprediction  of  the  peak 
streamwise  velocity  fluctuations.  The  peak  value  decreases  slightly  as  the  grid  is  refined. 
This  is  also  observed  for  the  DMM  and  the  DRM  (N=5).  It  should  be  noted  that  the 
DRM  actually  predicts  a peak  value  that  is  lower  than  the  DNS  data.  This  is  very  unusual 
in  LES,  because  most  models  will  overpredict  the  streamwise  velocity  fluctuations  and 
underpredict  the  wall-normal  and  spanwise  fluctuations.  However,  the  wall-normal  and 
spanwise  velocity  fluctuations  shown  in  figure  4 are  even  further  underpredicted  when 
applying  DRM  (N=5)  compared  to  the  other  two  models. 

The  modeled  shear  stresses  are  shown  in  figure  5.  It  is  a well-known  problem  that  the 
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Figure  4.  Velocity  fluctuations  in  the  streamwise  In'l,  wall-normal  and  spanwise  |ii;'|  di- 
rections. o;  filtered  DNS, : DSM  (64,49,48), : DSM  (128,97,96), ; DMM 

(64,49,48), : DMM  (128,97,96), : DRM  (N=5)  (64,49,48),  and : DRM  (N=5) 

(128,97,96). 


Figure  5.  Modeled  shear  stress,  ri2,  using  different  turbulence-closure  models,  o;  filtered  DNS, 

: DSM  (64,49,48), : DSM  (128,97,96), : DMM  (64,49,48), ; DMM 

(128,97,96), : DRM  (N=5)  (64,49,48),  and ; DRM  (N=5)  (128,97,96). 


DSM  does  not  predict  sufficient  shear  stress  in  the  near- wall  region  (Baggett,  Jimenez 
& Kravchenko  1997).  As  shown  in  the  figure,  the  DSM  predicts  the  lowest  peak  values, 
while  the  largest  are  produced  by  the  DRM  (N=5).  The  modeled  shear  stress  increases 
when  the  level  of  reconstruction  is  increased.  However,  the  contribution  from  the  DSM 
does  not  change  much  between  the  different  simulations;  the  peak  value  is  approximately 
the  same.  The  increase  of  modeled  shear  stress  is  therefore  almost  entirely  due  to  the 
RSFS  model.  The  peak  value  of  the  modeled  shear  stress  is  approaching  the  filtered 
DNS  data  as  the  level  of  reconstruction  is  increased.  Towards  the  center  of  the  channel, 
the  high  level  reconstruction  model  overpredicts  the  shear  stress  when  compared  to  the 
filtered  DNS  data. 
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The  true  LES  approach  is  investigated  in  turbulent  channel  flow  using  commutative  fil- 
ter functions  in  all  three  spatial  directions.  In  the  true  LES  approach,  a grid-independent 
solution  to  the  filtered  governing  equations  is  obtained.  The  LES  solution  depends  upon 
the  explicit  filter  width  used,  but  is  independent  of  the  computational  grid.  A compu- 
tational code  using  an  energy-conserving  fourth-order  finite-difference  scheme  is  applied 
and  fourth-order  commutative  filters  are  used.  Simulations  of  turbulent  channel  flow  were 
performed  at  Rer  = 395.  The  explicit  filter  width  was  kept  fixed  while  the  computational 
grid  was  refined,  to  obtain  a grid-independent  solution.  The  results  using  two  different 
grid  resolutions  show  only  minor  differences,  indicating  that  the  LES  solutions  are  nearly 
converged.  The  explicit  filtering  also  reduces  the  numerical  errors  that  are  associated 
with  the  high-wavenumber  portion  of  the  spectrum  when  using  non-spectral  methods. 
Therefore,  explicit  filtering  in  LES,  using  high-order  commutative  filters,  results  in  a 
numerically-clean  environment  where  turbulence-closure  models  can  be  investigated  in 
grid-independent  LES  solutions.  This  could  not  be  performed  using  the  traditional  LES 
approach,  since  the  contribution  from  the  SFS  models  decreases  as  the  computational 
grid  is  refined. 

The  turbulence-closure  models  investigated  are  the  DSM,  DMM  and  DRM.  The  mod- 
els are  compared  to  filtered  DNS  data  for  mean  velocity  profiles,  velocity  fluctuations, 
and  modeled  shear  stresses.  The  mean  velocity  profiles  and  the  streamwise  velocity  fluc- 
tuations improve  as  the  level  of  reconstruction  increases.  The  closest  agreement  between 
the  LES  results  and  the  filtered  DNS  data,  in  this  study,  is  obtained  by  the  level  of 
reconstruction  of  five  (DRM,  N=5). 

The  poor  agreement  between  the  filtered  DNS  results  and  the  DSM  shows  the  need 
for  RSFS  models  when  using  a smooth  (in  spectral  space)  explicit  filter  function.  The  re- 
sults predicted  by  the  models  investigated  show  a distinct  improvement  in  the  predicted 
quantities,  when  compared  to  filtered  DNS  results,  as  the  level  of  reconstruction  is  in- 
creased. These  improvements  are  probably  due  to  the  increase  of  modeled  shear  stress 
in  the  near-wall  region.  The  DSM  is  known  not  to  predict  sufiicient  shear  stress  in  the 
near-wall  region,  and  as  the  level  of  reconstruction  is  increased  so  is  the  modeled  shear 
stress.  The  increase  is  almost  entirely  due  to  the  RSFS  model,  since  the  contribution 
from  the  DSM  does  not  change  significantly. 
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